#
#    Copyright (C) 2016 Francisco Javier Parra <franparpe@openmailbox.org>
#
#    This program is free software; you can redistribute it and/or modify
#    it under the terms of the GNU General Public License as published by
#    the Free Software Foundation; either version 2 of the License, or
#    (at your option) any later version.
#
#    This program is distributed in the hope that it will be useful,
#    but WITHOUT ANY WARRANTY; without even the implied warranty of
#    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#    GNU General Public License for more details.
#
#    You should have received a copy of the GNU General Public License
#    along with this program; if not, write to the Free Software
#    Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301 USA
#

import math
import numpy
from qtcam.core.checkOutline import *


class SectionParameters(object):  # Because every cam uses the same thing
    """The paramaters that are present in every cam saction

    """
    def __init__(self, LowestPosition,
                 HighestPosition, StartAngle,
                 EndAngle, Precision, AngularVelocity):
        self.LowestPosition = LowestPosition
        self.HighestPosition = HighestPosition
        self.StartAngle = StartAngle
        self.EndAngle = EndAngle
        self.Precision = Precision
        self.AngularVelocity = AngularVelocity
        self.abcissaeValues = numpy.arange(self.StartAngle, self.EndAngle,
                                           self.Precision)
        # Beta is in degrees
        self.beta = self.StartAngle - self.EndAngle

    def ReCalculateAbcissaeValues(self):
        """Initialize the abcissaeValues class variable
(Better used in constant width cam)"""
        self.abcissaeValues = numpy.arange(self.StartAngle,
                                           self.EndAngle, self.Precision)
        
    def CreateYaxisPosition(self,):  # eliminar para crear algo mejor
        """This method initializes two arrays atributes, the first is the
        angle and the second is the Position"""
        self.ordinatePosition = self.ArrayPosition()

    def CreateYaxisVelocity(self,):
        """This method creates the ordinate values for the Velocity"""
        self.ordinateVelocity = self.ArrayVelocity()

    def CreateYaxisAcceleration(self,):
        """This method creates the ordinate values for the Acceleration"""
        self.ordinateAcceleration = self.ArrayAcceleration()

    def CreateYaxisOverAcceleration(self,):
        """This method creates the ordinate values for the
        overacceleration"""
        self.ordinateOverAcceleration = self.ArrayOverAcceleration()

    def generate_cam(self, rb=0, excentricity=0):
        "This method generates the points that makes the cam"
        self.CreateYaxisPosition()
        pointList = []
        # This loop creates all points that generates the cam
        for angle, postion in zip(self.abcissaeValues, self.ordinatePosition):
            # "S" is the matrix that changes the basis and will be
            # divided in two columns: SC1 and SC2
            cosFI = math.cos(math.radians(angle))
            sinFI = math.sin(math.radians(angle))
            x = (cosFI * excentricity) + (sinFI * (postion + rb))
            y = ((-sinFI) * excentricity) + (cosFI * (postion + rb))
            point = (x, y)
            pointList.append(point)

        return pointList

    def generate_cam_flat_follower(self, rb, excentricity=0, beta=0):
        """Generates all points for the polygon to be constructed.

        rb: the base radius of the cam.
        
        excentricity: the perpendicular distance the follower
        has with the center of rotation of the cam.

        beta: the inclination in degrees that the follower have with the vertical.

        Keyword Arguments:
        rb           --
        excentricity -- (default 0)
        beta         -- (default 0)
        """
        curvatureRadiusList = []
        # This list will contain all angles whose radius value is zero
        anglesE0Radius = []
        pointsE0Radius = []
        # This list will contain all angles whose radius value
        # is less than zero
        anglesL0Radius = []
        pointsL0Radius = []
        cosBeta = math.cos(math.radians(beta))
        sinBeta = math.sin(math.radians(beta))
        initDistance = (rb + (excentricity * sinBeta)) / cosBeta
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        self.CreateYaxisAcceleration()
        pointList = []
        for angle, position, deriv, deriv2 in zip(self.abcissaeValues,
                                                  self.ordinatePosition,
                                                  self.ordinateVelocity,
                                                  self.ordinateAcceleration):
            deriv = deriv / self.AngularVelocity
            deriv2 = deriv2 / math.pow(self.AngularVelocity, 2)
            cosFI = math.cos(math.radians(angle))
            sinFI = math.sin(math.radians(angle))
            # aplying change of basis
            x = (((deriv * cosBeta) * cosFI) +
                 ((((position + initDistance) * cosBeta) -
                   (excentricity * sinBeta)) * sinFI))
            y = ((((deriv * cosBeta) * (-sinFI)) +
                 ((((position + initDistance) * cosBeta) -
                   (excentricity * sinBeta)) * cosFI)))
            point = (x, y)
            pointList.append(point)
            # checking if the outline is OK by checking the curvature radius
            curvatureRadius = ((position + initDistance + deriv2) *
                               cosBeta) - (excentricity * sinBeta)
            if (curvatureRadius == 0):  # the outline has a peak
                anglesE0Radius.append(curvatureRadius)
                pointsE0Radius.append(point)
            elif (curvatureRadius < 0):
                anglesL0Radius.append(curvatureRadius)
                pointsL0Radius.append(point)
            else:
                pass
            curvatureRadiusList.append(curvatureRadius)
        # Return the possible faliures
        # ShowCurvatureRadiusFailures(anglesE0Radius,
        #                             anglesL0Radius,
        #                             pointsE0Radius,
        #                             pointsL0Radius)
        return pointList, curvatureRadiusList

    def generate_cam_circular_follower(self, rb,
                                       followerRadius, excentricity=0):
        """Generates the points that generates the outline"""
        curvatureRadiusPList = []
        curvatureRadiusList = []
        pressureAngleList = []
        offset = []
        # This list will contain all angles correspondeing to a peak
        peakAngles = []
        peakPoints = []
        # This list will contain all angles where there is no contact
        # with the follower
        noContactAngles = []
        noContactPoints = []
        # Angulos con rebaje
        negAngles = []
        negPoints = []
        # Angles with no good pressure angle
        badPressureAngles = []
        badPressurePoints = []
        pressureAngles = []
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        self.CreateYaxisAcceleration()
        initDistance = math.sqrt(math.pow((rb + followerRadius), 2) -
                                 math.pow(excentricity, 2))
        pointList = []
        for angle, position, velocity, acceleration in zip(self.abcissaeValues,
                                                           self.ordinatePosition,
                                                           self.ordinateVelocity,
                                                           self.ordinateAcceleration):
            velocity = velocity / self.AngularVelocity
            acceleration = acceleration / math.pow(self.AngularVelocity, 2)
            # Whe have a rotating basis so we multiply the two components
            # of the vector with the change of basis matrix
            # [S] = [[cosFI, sinFI];[-sinFI, cosFI]]
            # the rotation is clockwise

            # d is the total distance from the center of the cam
            # to the center of the follower's wheel
            d = position + initDistance
            cosFI = math.cos(math.radians(angle))
            sinFI = math.sin(math.radians(angle))
            OCx = (cosFI * excentricity) + (sinFI * d)
            OCy = ((-sinFI) * excentricity) + (cosFI *
                                               (position + initDistance))
            offset.append((OCx, OCy))
            # OC is the point that goes from the center of the cam to the
            # center of the follower's wheel
            # We calculate tangent vector
            vectorTanx = (((velocity - excentricity) * sinFI) +
                          ((position + initDistance) * cosFI))
            vectorTany = (((velocity - excentricity) * cosFI) -
                          ((position + initDistance) * sinFI))
            # We now have to rotate the tangent vector -90 degrees to obtain
            # the normal vector. The rotation matrix for a counterclockwise
            # spin is: [S] = [[0, 1];[-1, 0]]
            vectorNormalx = vectorTany
            vectorNormaly = -vectorTanx
            # We have interest in the unit vector
            vectorNormalUnitx = (vectorNormalx /
                                 math.hypot(vectorNormalx, vectorNormaly))
            vectorNormalUnity = (vectorNormaly /
                                 math.hypot(vectorNormalx, vectorNormaly))
            # we now calculate OJ, the vector that points
            # to every point that generates the outline
            OJx = (OCx + followerRadius * vectorNormalUnitx)
            OJy = (OCy + followerRadius * vectorNormalUnity)
            OJ = (OJx, OJy)
            pointList.append(OJ)

            # Here comes the radius of curvature calculation
            curvatureRadiusP = (math.pow((math.pow(d, 2) +
                                          math.pow((velocity - excentricity), 2)), 3/2) /
                                (((velocity - excentricity) *
                                 (2 * velocity - excentricity)) -
                                 (d *
                                  (acceleration - (d)))))
            curvatureRadius = curvatureRadiusP - followerRadius
            if (curvatureRadius < 0):
                if not (abs(curvatureRadius) > followerRadius):
                    noContactAngles.append(angle)
                    noContactPoints.append(OJ)
                else:
                    pass
            elif (curvatureRadiusP > 0):
                if (curvatureRadiusP == followerRadius):  # There is a peak
                    peakAngles.append(angle)
                    noContactPoints.append(OJ)
                elif (curvatureRadiusP < followerRadius):  # Rebaje
                    negAngles.append(angle)
                    negPoints.append(OJ)
                else:
                    pass
                
            curvatureRadiusPList.append(curvatureRadiusP)
            curvatureRadiusList.append(curvatureRadius)

            # Checking by pressure angle
            pAngle = math.atan((velocity - excentricity) /
                               (d + math.sqrt(math.pow((rb + followerRadius), 2) +
                                              math.pow(excentricity, 2))))  # Radiants
            if abs(math.degrees(pAngle)) > 30:
                pressureAngles.append(pAngle)
                badPressureAngles.append(math.degrees(angle))
                badPressurePoints.append(OJ)
            pressureAngleList.append(math.degrees(pAngle))

        # ShowPressureAngleFailuresCircular(pressureAngles,
        #                                   badPressureAngles, badPressurePoints)
        # ShowCurvatureRadiusFailuresCircular(peakAngles, negAngles,
        #                                     noContactAngles,
        #                                     peakPoints,
        #                                     negPoints,
        #                                     noContactPoints)

        return pointList, curvatureRadiusPList, curvatureRadiusList, pressureAngleList, offset

    def generate_cam_rotation_follower(self, rb, l1, l2):
        """
        Keyword Arguments:
        l1 -- 
        l2 -- 
        """
        numerator = math.pow(l1, 2) + math.pow(l2, 2) - math.pow(rb, 2)
        denominator = 2 * l1 * l2
        # Fi2o is the initial angle the follower has because of the base
        # radius of the cam
        Fi2o = math.degrees(math.acos(numerator / denominator))
        print(Fi2o)  # Trying
        self.CreateYaxisPosition()
        outline = []
        for Fi1, ordinateAngle in zip(self.abcissaeValues,
                                      self.ordinatePosition):
            Fi2 = ordinateAngle + Fi2o
            # print(Fi2)  # Trying
            sinFi2 = math.sin(math.radians(Fi2))
            cosFi2 = math.cos(math.radians(Fi2))
            x = l1 - (l2 * (cosFi2))
            y = l2 * sinFi2
            # OJ is the vector that has its origin
            # in the center of the cam and is pointing at the
            # perimeter of the cam with the fixed base
            # To change the basis to the fixed one, We have to multiply
            # the current vector with the the following matrix:
            # [S] = [[cosFi1, sinFi1];[-sinFi1, cosFi1]]
            # Operation: [S][OJ]rotateBasis
            sinFi1 = math.sin(math.radians(Fi1))
            cosFi1 = math.cos(math.radians(Fi1))
            OJx = ((x * cosFi1) + (y * sinFi1))
            OJy = ((x * (-sinFi1)) + (y * cosFi1))
            OJ = (OJx, OJy)
            outline.append(OJ)

        return outline

    def generate_cam_rotation_follower_plane(self, rb, l1, l2):
        """Generates the outline for the cam with an angular follower
        with a plane instead of a point as the contact part
        Keyword Arguments:
        rb -- 
        l1 -- 
        l2 --
        """
        curvatureRadiusList = []
        # This list will contain all angles whose radius value is zero
        anglesE0Radius = []
        pointsE0Radius = []
        # This list will contain all angles whose radius value
        # is less than zero
        anglesL0Radius = []
        pointsL0Radius = []

        LowestPositionDegrees = self.LowestPosition
        HighestPositionDegrees = self.HighestPosition
        self.LowestPosition = math.radians(self.LowestPosition)
        self.HighestPosition = math.radians(self.HighestPosition)
        # Passing values to radiants
        outline = []
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        self.CreateYaxisAcceleration()
        Fi2o = math.asin((rb - l2) / l1)
        for Fi1, ordinateAngle, ordinateDeriv, ordinateDeriv2 in zip(self.abcissaeValues,
                                                                     self.ordinatePosition,
                                                                     self.ordinateVelocity,
                                                                     self.ordinateAcceleration):
            ordinateDeriv = ordinateDeriv / self.AngularVelocity
            ordinateDeriv2 = ordinateDeriv2 / math.pow(self.AngularVelocity, 2)
            Fi2 = Fi2o + ordinateAngle
            sinFi2 = math.sin(Fi2)
            cosFi2 = math.cos(Fi2)
            l3 = (l1 * cosFi2) / (1 + ordinateDeriv)
            OJ1 = (l1 + (l2 * sinFi2) - (l3 * cosFi2))
            OJ2 = ((l2 * cosFi2) + (l3 * sinFi2))
            # Time to change the basis to the fix one
            sinFi1 = math.sin(math.radians(Fi1))
            cosFi1 = math.cos(math.radians(Fi1))
            OJx = ((OJ1 * cosFi1) + (OJ2 * sinFi1))
            OJy = ((OJ1 * (-sinFi1)) + (OJ2 * cosFi1))
            OJ = (OJx, OJy)
            outline.append(OJ)

            # For this cam, only is necessary to confirm
            # if the Radius of curvature is correct
            firstAddend = ((1 + (2 * ordinateDeriv)) /
                           math.pow((1 + ordinateDeriv), 2)) * sinFi2
            secondAddend = (ordinateDeriv2 /
                            math.pow((1 + ordinateDeriv), 3)) * cosFi2
            giantExpression = firstAddend + secondAddend
            radiusOfCurvature = l2 + (l1 * giantExpression)
            if (radiusOfCurvature == 0):  # the outline has a peak
                anglesE0Radius.append(radiusOfCurvature)
                pointsE0Radius.append(OJ)
            elif (radiusOfCurvature < 0):
                anglesL0Radius.append(radiusOfCurvature)
                pointsL0Radius.append(OJ)
            else:
                pass
            curvatureRadiusList.append(radiusOfCurvature)
        # Return the possible faliures
        # ShowCurvatureRadiusFailures(anglesE0Radius,
        #                             anglesL0Radius,
        #                             pointsE0Radius,
        #                             pointsL0Radius)
            
        self.LowestPosition = LowestPositionDegrees
        self.HighestPosition = HighestPositionDegrees
        return outline, curvatureRadiusList

    def generate_cam_rotation_follower_circular(self, rb, rf, l1, l2):
        """
        arguments:
        rb: The base radius of the cam,

        rf: The radius of the circular follower

        l1: The distance...

        l2: The longitude of the arm.
        """
        # The lists
        curvatureRadiusPList = []
        curvatureRadiusList = []
        pressureAngleList = []
        # keeping things to genereate
        LowestPositionDegrees = self.LowestPosition
        HighestPositionDegrees = self.HighestPosition
        self.LowestPosition = math.radians(self.LowestPosition)
        self.HighestPosition = math.radians(self.HighestPosition)
        # This list will contain all angles correspondeing to a peak
        peakAngles = []
        peakPoints = []
        # This list will contain all angles where there is no contact
        # with the follower
        noContactAngles = []
        noContactPoints = []
        # Angulos con rebaje
        negAngles = []
        negPoints = []
        # Angles with no good pressure angle
        badPressureAngles = []
        badPressurePoints = []
        # All pressure angles that are problematic
        pressureAngles = []
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        self.CreateYaxisAcceleration()
        outline = []
        offset = []
        # Fi1 is the angle that symbolize the rotation of the cam
        # Fi2 is the angle the follower does with the horizontal, Fi2
        # depends on Fi1
        numerator = math.pow(l1, 2) + math.pow(l2, 2) - math.pow(rb + rf, 2)
        denominator = 2 * l1 * l2
        # Fi2o is the initial angle the follower has because of the base
        # radius of the cam
        Fi2o = math.acos(numerator / denominator)
        for Fi1, ordinateAngle, angleDerivative, angleDerivative2 in zip(self.abcissaeValues,
                                                                         self.ordinatePosition,
                                                                         self.ordinateVelocity,
                                                                         self.ordinateAcceleration):
            angleDerivative = angleDerivative / self.AngularVelocity
            angleDerivative2 = angleDerivative2 / math.pow(self.AngularVelocity, 2)
            Fi2 = Fi2o + ordinateAngle
            cosFi2 = math.cos(Fi2)
            sinFi2 = math.sin(Fi2)
            # Vector OC calculation (the one that starts on the center of
            # the cam and ends on the center of the circular follower).
            # The initial basis is moving with the rotation
            OC1 = l1 - (l2 * cosFi2)
            OC2 = l2 * sinFi2
            # Now we change the basis of the vector
            Fi1 = math.radians(Fi1)
            cosFi1 = math.cos(Fi1)
            sinFi1 = math.sin(Fi1)
            OCx = (cosFi1 * OC1) + (sinFi1 * OC2)
            OCy = ((-sinFi1) * OC1) + (cosFi1 * OC2)
            offset.append((OCx, OCy))
            
            # Tangentk vector calculation
            sinFi1PlusFi2 = math.sin(Fi1 + Fi2)
            cosFi1PlusFi2 = math.cos(Fi1 + Fi2)
            
            vectorTanx = (((-l1) * sinFi1) +
                          (l2 * sinFi1PlusFi2 * (1 + angleDerivative)))
            vectorTany = (((-l1) * cosFi1) +  # Había un error aquí en la publicación de Zayas
                          (l2 * cosFi1PlusFi2 * (1 + angleDerivative)))
            
            # We now have to rotate the tangent vector -90 degrees to obtain
            # the normal vector. The rotation matrix for a counterclockwise
            # spin is: [S] = [[0, 1];[-1, 0]]
            vectorNormalx = vectorTany
            vectorNormaly = -vectorTanx
            vectorNormalUnitx = (vectorNormalx /
                                 math.hypot(vectorNormalx, vectorNormaly))
            vectorNormalUnity = (vectorNormaly /
                                 math.hypot(vectorNormalx, vectorNormaly))
            # we now calculate OJ, the vector that points
            # to every point that generates the outline
            OJx = (OCx + (rf * vectorNormalUnitx))
            OJy = (OCy + (rf * vectorNormalUnity))
            OJ = (OJx, OJy)  # testing
            outline.append(OJ)  # testing

            # Cheking the radius of curvature
            C = (-l1 + (l2 * math.cos(Fi2) * (1 + angleDerivative)))  # here
            D = l2 * math.sin(Fi2) * (1 + angleDerivative)
            Cp2 = math.pow(C, 2)
            Dp2 = math.pow(D, 2)
            curvatureRadiusP = (math.pow((Cp2 + Dp2), (3 / 2)) /  # numerator
                                ((((C * l2 *
                                    angleDerivative * cosFi2) + Dp2) *
                                  (1 + angleDerivative)) +  # +
                                 Cp2 -  # -
                                 (l1 * sinFi2 * l2 * angleDerivative2)))
            # curvatureRadiusP = (math.pow((Cp2 + Dp2), (3 / 2)) /
            #                     (((Cp2 + Dp2) * (1 + angleDerivative)) -
            #                      (l1 * C * angleDerivative) +
            #                     (l1 * l2 * sinFi2 * angleDerivative2)))
            # curvatureRadiusP = (math.pow((Cp2 + Dp2), (3 / 2)) /
            #                     ((C * l2 * angleDerivative2 * sinFi2) +
            #                     (C * l2 *
            #                      (angleDerivative + 1) *
            #                      angleDerivative * cosFi2) +
            #                      Cp2 -
            #                      (D * l2 * angleDerivative2 * cosFi2) +
            #                      (angleDerivative * Dp2) +
            #                      Dp2))
            curvatureRadius = curvatureRadiusP - rf
            #print(curvatureRadius)
            #print("P--->", curvatureRadiusP)
            #print(C)
            #print(D)
            #print(math.pow((Cp2 + Dp2), (3 / 2)))
            if (curvatureRadius < 0):
                if abs(curvatureRadius) <= rf:
                    noContactAngles.append(math.degrees(Fi1))
                    noContactPoints.append(OJ)
                else:
                    pass
            elif (curvatureRadiusP > 0):
                if (curvatureRadiusP == rf):  # There is a peak
                    peakAngles.append(math.degrees(Fi1))
                    peakPoints.append(OJ)
                elif (curvatureRadiusP < rf):  # Rebaje
                    negAngles.append(math.degrees(Fi1))
                    negPoints.append(OJ)
                else:
                    pass
            curvatureRadiusPList.append(curvatureRadiusP)
            curvatureRadiusList.append(curvatureRadius)
            # Time to calculate the pressure angle
            OCmodulus = math.hypot(OC1, OC2)
            numerator = (math.pow(OCmodulus, 2) +
                         math.pow(l2, 2) -
                         math.pow(l1, 2))
            denominator = 2 * OCmodulus * l2
            gamma = math.acos(numerator / denominator)
            # numerator = (l2 * sinFi2 * angleDerivative)
            # denominator = (OCmodulus - numerator)
            numerator = (l2 * math.sin(gamma) * angleDerivative)
            denominator = (OCmodulus -
                           (l2 * math.cos(gamma) * angleDerivative))
            alpha = math.atan(numerator / denominator)
            pressureAngle = gamma + alpha - (math.pi / 2)
            if (abs(math.degrees(pressureAngle)) > 30):
                pressureAngles.append(pressureAngle)
                badPressureAngles.append(math.degrees(Fi1))
                badPressurePoints.append(OJ)
            else:
                pass
            pressureAngleList.append(math.degrees(pressureAngle))
        # Falta el tratamiento final
        # ShowCurvatureRadiusFailuresCircular(peakAngles,
        #                                     negAngles,
        #                                     noContactAngles,
        #                                     peakPoints,
        #                                     negPoints,
        #                                     noContactPoints)
        # ShowPressureAngleFailuresCircular(pressureAngles,
        #                                   badPressureAngles, badPressurePoints)
        self.LowestPosition = LowestPositionDegrees
        self.HighestPosition = HighestPositionDegrees

        return (outline,
                curvatureRadiusPList,
                curvatureRadiusList,
                pressureAngleList,
                offset)

    def Conjugate_plain(self, rb, beta, excentricity, separation):
        """Returns the outline as a list of points of a demodromic cam
        
        arguments:
        separation: A float that represents the distance between followers

        excentricity: A float that represents the excentricity
        """
        # all curvRadius
        curvatureRadiusListUp = []
        curvatureRadiusListLo = []
        # This list will contain all angles whose radius value is zero
        radiusOfCurvE0Sup = []
        pointsE0RadiusSup = []
        # This list will contain all angles whose radius value
        # is less than zero
        radiusOfCurvL0Sup = []
        pointsL0RadiusSup = []
        
        # This list will contain all angles whose radius value is zero
        radiusOfCurvE0Inf = []
        pointsE0RadiusInf = []
        # This list will contain all angles whose radius value
        # is less than zero
        radiusOfCurvL0Inf = []
        pointsL0RadiusInf = []
        cosBeta = math.cos(math.radians(beta))
        sinBeta = math.sin(math.radians(beta))
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        self.CreateYaxisAcceleration()
        sinBeta = math.sin(math.radians(beta))
        cosBeta = math.cos(math.radians(beta))
        initDistance = (rb + (excentricity * sinBeta)) / cosBeta
        outline1 = []
        outline2 = []
        for fi1, position, deriv, deriv2 in zip(self.abcissaeValues,
                                                self.ordinatePosition,
                                                self.ordinateVelocity,
                                                self.ordinateAcceleration):
            deriv = deriv / self.AngularVelocity
            deriv2 = deriv2 / math.pow(self.AngularVelocity, 2)
            cosFI = math.cos(math.radians(fi1))
            sinFI = math.sin(math.radians(fi1))
            d1 = initDistance + position
            d2 = d1 - separation
            # aplying change of basis
            x = (((deriv * cosBeta) * cosFI) +
                 (((d1 * cosBeta) -
                   (excentricity * sinBeta)) * sinFI))
            y = ((((deriv * cosBeta) * (-sinFI)) +
                 (((d1 * cosBeta) -
                   (excentricity * sinBeta)) * cosFI)))
            point1 = (x, y)
            x2 = (((deriv * cosBeta) * cosFI) +
                  (((d2 * cosBeta) -
                   (excentricity * sinBeta)) * sinFI))
            y2 = ((((deriv * cosBeta) * (-sinFI)) +
                   (((d2 * cosBeta) -
                     (excentricity * sinBeta)) * cosFI)))
            point2 = (x2, y2)
            outline2.append(point2)
            outline1.append(point1)
            # checking if the outline is OK by checking the curvature radius
            radiusOfCurvatureSup = ((d1 + deriv2) *
                                    cosBeta) - (excentricity * sinBeta)
            # Less sign is requred because d2 is in the other side
            radiusOfCurvatureInf = ((-d2 + deriv2) *  # Less sign is requred
                                    cosBeta) - (excentricity * sinBeta)
            if (radiusOfCurvatureSup == 0):  # the outline has a peak
                radiusOfCurvE0Sup.append(radiusOfCurvatureSup)
                pointsE0RadiusSup.append(point1)
            elif (radiusOfCurvatureSup < 0):
                radiusOfCurvL0Sup.append(radiusOfCurvatureSup)
                pointsL0RadiusSup.append(point1)
            else:
                pass

            if (radiusOfCurvatureInf == 0):  # the outline has a peak
                radiusOfCurvE0Inf.append(radiusOfCurvatureInf)
                pointsE0RadiusInf.append(point2)
            elif (radiusOfCurvatureInf < 0):
                radiusOfCurvL0Inf.append(radiusOfCurvatureInf)
                pointsL0RadiusInf.append(point2)
            else:
                pass
            curvatureRadiusListUp.append(radiusOfCurvatureSup)
            curvatureRadiusListLo.append(radiusOfCurvatureInf)
        # Return the Failures
        ShowCurvatureRadiusFailuresConjugatePlane(radiusOfCurvE0Sup,
                                                  radiusOfCurvL0Sup,
                                                  pointsE0RadiusSup,
                                                  pointsL0RadiusSup,
                                                  radiusOfCurvE0Inf,
                                                  radiusOfCurvL0Inf,
                                                  pointsE0RadiusInf,
                                                  pointsL0RadiusInf)
        return outline1, outline2, curvatureRadiusListUp, curvatureRadiusListLo

    def Conjugate_circular(self, rb, followerRadius, excentricity, separation):
        """Documentation goes here"""
        # curvature radius list
        curvatureRadiusPListUp = []
        curvatureRadiusPListLo = []
        curvatureRadiusListUp = []
        curvatureRadiusListLo = []

        offset1 = []
        offset2 = []
        # pressure angle list
        pressureAngleListUp = []
        pressureAngleListLo = []
        # This list will contain all angles correspondeing to a peak
        peakAngles = []
        peakPoints = []
        # This list will contain all angles where there is no contact
        # with the follower
        noContactAngles = []
        noContactPoints = []
        # Angulos con rebaje
        negAngles = []
        negPoints = []
        # Angles with no good pressure angle
        pressureAngles = []
        badPressureAngles = []
        badPressurePoints = []
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        self.CreateYaxisAcceleration()
        initDistance = math.sqrt(math.pow((rb + followerRadius), 2) -
                                 math.pow(excentricity, 2))
        outline1 = []
        outline2 = []
        for fi1, position, deriv, deriv2 in zip(self.abcissaeValues,
                                                self.ordinatePosition,
                                                self.ordinateVelocity,
                                                self.ordinateAcceleration):
            deriv = deriv / self.AngularVelocity
            deriv2 = deriv2 / math.pow(self.AngularVelocity, 2)
            # Whe have a rotating basis so we multiply the two components
            # of the vector with the change of basis matrix
            # [S] = [[cosFI, sinFI];[-sinFI, cosFI]]
            # the rotation is clockwise
            d = position + initDistance
            d2 = d - separation
            cosFI = math.cos(math.radians(fi1))
            sinFI = math.sin(math.radians(fi1))
            OCx1 = (cosFI * excentricity) + (sinFI * d)
            OCy1 = ((-sinFI) * excentricity) + (cosFI * d)

            OCx2 = ((cosFI * excentricity) +
                    (sinFI * (position + (initDistance - separation))))
            OCy2 = (((-sinFI) * excentricity) +
                    (cosFI * (position + (initDistance - separation))))
            offset1.append((OCx1, OCy1))
            offset2.append((OCx2, OCy2))
            # OC is the point that goes from the center of the cam to the
            # center of the follower's wheel
            # We calculate tangent vector
            vectorTanx = (((deriv - excentricity) * sinFI) +
                          ((position + initDistance) * cosFI))
            vectorTany = (((deriv - excentricity) * cosFI) -
                          ((position + initDistance) * sinFI))
            vectorTanx2 = (((deriv - excentricity) * sinFI) +
                           ((position + initDistance - separation) * cosFI))
            vectorTany2 = (((deriv - excentricity) * cosFI) -
                           ((position + initDistance - separation) * sinFI))
            # We now have to rotate the tangent vector -90 degrees to obtain
            # the normal vector. The rotation matrix for a counterclockwise
            # spin is: [S] = [[0, 1];[-1, 0]]
            vectorNormalx1 = vectorTany
            vectorNormaly1 = -vectorTanx
            vectorNormalx2 = vectorTany2
            vectorNormaly2 = -vectorTanx2
            # We have interest in the unit vector
            vectorNormalUnitx1 = (vectorNormalx1 /
                                  math.hypot(vectorNormalx1, vectorNormaly1))
            vectorNormalUnity1 = (vectorNormaly1 /
                                  math.hypot(vectorNormalx1, vectorNormaly1))
            vectorNormalUnitx2 = (vectorNormalx2 /
                                  math.hypot(vectorNormalx2, vectorNormaly2))
            vectorNormalUnity2 = (vectorNormaly2 /
                                  math.hypot(vectorNormalx2, vectorNormaly2))
            # we now calculate OJ, the vector that points
            # to every point that generates the outline
            OJx1 = (OCx1 + followerRadius * vectorNormalUnitx1)
            OJy1 = (OCy1 + followerRadius * vectorNormalUnity1)
            OJx2 = (OCx2 + followerRadius * vectorNormalUnitx2)
            OJy2 = (OCy2 + followerRadius * vectorNormalUnity2)
            OJ1 = (OJx1, OJy1)
            OJ2 = (OJx2, OJy2)
            outline1.append(OJ1)
            outline2.append(OJ2)
            
            # Here comes the radius of curvature calculation
            # The upper part? Is it necessary?
            curvatureRadiusP = (math.pow((math.pow(d, 2) +
                                          math.pow((deriv - excentricity), 2)), 3/2) /
                                ((deriv - excentricity) * (2 * deriv - excentricity) -
                                 (d *
                                  (deriv2 - d))))
            curvatureRadius = curvatureRadiusP - followerRadius
            if (curvatureRadius < 0):
                if not (abs(curvatureRadius) > followerRadius):
                    noContactAngles.append(angle)
                    noContactPoints.append((OJ1, OJ2))
                else:
                    pass
            elif (curvatureRadiusP > 0):
                if (curvatureRadiusP == followerRadius):  # There is a peak
                    peakAngles.append(angle)
                    noContactPoints.append((OJ1, OJ2))
                elif (curvatureRadiusP < followerRadius):  # Rebaje
                    negAngles.append(angle)
                    negPoints.append((OJ1, OJ2))
                else:
                    pass
            curvatureRadiusPListUp.append(curvatureRadiusP)
            curvatureRadiusListUp.append(curvatureRadius)


            # Here comes the radius of curvature calculation
            # The lower part? Is it necessary?
            curvatureRadiusPInf = (math.pow((math.pow(d2, 2) +
                                             math.pow((deriv - excentricity),
                                                      2)), 3/2) /
                                   ((deriv - excentricity) *
                                    (2 * deriv - excentricity) -
                                    (d2 *
                                     (deriv2 - d2))))
            curvatureRadiusInf = curvatureRadiusPInf - followerRadius
            if (curvatureRadius < 0):
                if not (abs(curvatureRadius) > followerRadius):
                    noContactAngles.append(angle)
                    noContactPoints.append((OJ1, OJ2))
                else:
                    pass
            elif (curvatureRadiusP > 0):
                if (curvatureRadiusP == followerRadius):  # There is a peak
                    peakAngles.append(angle)
                    noContactPoints.append((OJ1, OJ2))
                elif (curvatureRadiusP < followerRadius):  # Rebaje
                    negAngles.append(angle)
                    negPoints.append((OJ1, OJ2))
                else:
                    pass
            curvatureRadiusPListLo.append(curvatureRadiusPInf)
            curvatureRadiusListLo.append(curvatureRadiusInf)
        

            # Checking by pressure angle upper
            pAngle = math.atan((deriv - excentricity) /
                               (d + math.sqrt(math.pow((rb + followerRadius),
                                                       2) +
                                              math.pow(excentricity, 2))))  # Radiants
            if abs(math.degrees(pAngle)) > 30:
                badPressureAngles.append(angle)
                badPressurePoints.append((OJ1, OJ2))
                pressureAngles.append(pAngle)

            # Lower pressure angle
            pAngleInf = math.atan((deriv - excentricity) /
                                  (-d2 + math.sqrt(math.pow((rb + followerRadius),
                                                            2) +
                                                   math.pow(excentricity, 2))))
            if abs(math.degrees(pAngleInf)) > 30:
                badPressureAngles.append(angle)
                badPressurePoints.append((OJ1, OJ2))
                pressureAngles.append(pAngle)

            pressureAngleListUp.append(math.degrees(pAngle))
            pressureAngleListLo.append(math.degrees(pAngleInf))

        return outline1, outline2, curvatureRadiusPListUp, curvatureRadiusPListLo, curvatureRadiusListUp, curvatureRadiusListLo, pressureAngleListUp, pressureAngleListLo, offset1, offset2

    def Conjugate_angular_plane(self, rb, beta, l1, l2sup, l2inf):
        """Returns the points of a conjugate cam... Blabla
        Keyword Arguments:
        rb    -- 
        l1    -- 
        l2sup -- 
        l2inf -- 
        beta -- The angle between
        """
        # The curvature radius list
        curvatureRadiusListUp = []
        curvatureRadiusListLo = []
        # This list will contain all bad curvature radius
        radiusOfCurvE0Sup = []
        pointsE0RadiusSup = []
        radiusOfCurvE0Inf = []
        pointsE0RadiusInf = []
        # This list will contain all angles whose radius value
        # is less than zero
        radiusOfCurvL0Sup = []
        pointsL0RadiusSup = []
        radiusOfCurvL0Inf = []
        pointsL0RadiusInf = []

        LowestPositionDegrees = self.LowestPosition
        HighestPositionDegrees = self.HighestPosition
        self.LowestPosition = math.radians(self.LowestPosition)
        self.HighestPosition = math.radians(self.HighestPosition)
        beta = math.radians(beta)
        # Hablar con Zayas, estoy usando la beta de la apertura
        outline1 = []
        outline2 = []
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        self.CreateYaxisAcceleration()
        Fi2o = math.asin((rb - l2sup) / l1)
        for Fi1, ordinateAngle, ordinateDeriv, ordinateDeriv2 in zip(self.abcissaeValues,
                                                                     self.ordinatePosition,
                                                                     self.ordinateVelocity,
                                                                     self.ordinateAcceleration):
            ordinateDeriv = ordinateDeriv / self.AngularVelocity
            ordinateDeriv2 = ordinateDeriv2 / math.pow(self.AngularVelocity, 2)
            Fi2_1 = Fi2o + ordinateAngle
            Fi2_2 = Fi2o + (ordinateAngle - beta)
            sinFi2_1 = math.sin(Fi2_1)
            cosFi2_1 = math.cos(Fi2_1)
            sinFi2_2 = math.sin(Fi2_2)
            cosFi2_2 = math.cos(Fi2_2)
            l3sup = ((l1 * cosFi2_1) / (ordinateDeriv + 1))
            l3inf = ((l1 * cosFi2_2) / (ordinateDeriv + 1))
            OJ1_1 = (l1 + (l2sup * sinFi2_1) - (l3sup * cosFi2_1))
            OJ2_1 = ((l2sup * cosFi2_1) + (l3sup * sinFi2_1))
            OJ1_2 = (l1 - (l2inf * sinFi2_2) - (l3inf * cosFi2_2))
            OJ2_2 = ((-l2inf * cosFi2_2) + (l3inf * sinFi2_2))
            sinFi1 = math.sin(math.radians(Fi1))
            cosFi1 = math.cos(math.radians(Fi1))
            OJx1 = ((OJ1_1 * cosFi1) + (OJ2_1 * sinFi1))
            OJy1 = ((OJ1_1 * (-sinFi1)) + (OJ2_1 * cosFi1))
            OJx2 = ((OJ1_2 * cosFi1) + (OJ2_2 * sinFi1))
            OJy2 = ((OJ1_2 * (-sinFi1)) + (OJ2_2 * cosFi1))
            point1 = (OJx1, OJy1)
            point2 = (OJx2, OJy2)
            outline1.append(point1)
            outline2.append(point2)
            # For this cam, only is necessary to confirm
            # if the Radius of curvature is correct
            firstAddend = ((1 + (2 * ordinateDeriv)) /
                           math.pow((1 + ordinateDeriv), 2)) * sinFi2_1
            secondAddend = (ordinateDeriv2 /
                            math.pow((1 + ordinateDeriv), 3)) * cosFi2_1
            giantExpression = firstAddend + secondAddend
            radiusOfCurvatureSup = l2sup + (l1 * giantExpression)

            curvatureRadiusListUp.append(radiusOfCurvatureSup)
            # Now for the inferior follower
            sinFi2_2 = math.sin(-Fi2_2)
            cosFi2_2 = math.cos(-Fi2_2)
            firstAddendInf = ((1 + (2 * ordinateDeriv)) /
                              math.pow((1 + ordinateDeriv), 2)) * sinFi2_2
            secondAddendInf = (ordinateDeriv2 /
                               math.pow((1 + ordinateDeriv), 3)) * cosFi2_2
            giantExpressionInf = firstAddendInf + secondAddendInf
            radiusOfCurvatureInf = l2inf + (l1 * giantExpressionInf)
            if (radiusOfCurvatureSup == 0):  # the outline has a peak
                radiusOfCurvE0Sup.append(radiusOfCurvatureSup)
                pointsE0RadiusSup.append(point1)
            elif (radiusOfCurvatureSup < 0):
                radiusOfCurvL0Sup.append(radiusOfCurvatureSup)
                pointsL0RadiusSup.append(point1)
            else:
                pass

            if (radiusOfCurvatureInf == 0):  # the outline has a peak
                radiusOfCurvE0Inf.append(radiusOfCurvatureInf)
                pointsE0RadiusInf.append(point2)
            elif (radiusOfCurvatureInf < 0):
                radiusOfCurvL0Inf.append(radiusOfCurvatureInf)
                pointsL0RadiusInf.append(point2)
            else:
                pass
            curvatureRadiusListLo.append(radiusOfCurvatureInf)

        # Return the possible faliures
        ShowCurvatureRadiusFailuresConjugatePlane(radiusOfCurvE0Sup,
                                                  radiusOfCurvL0Sup,
                                                  pointsE0RadiusSup,
                                                  pointsL0RadiusSup,
                                                  radiusOfCurvE0Inf,
                                                  radiusOfCurvL0Inf,
                                                  pointsE0RadiusInf,
                                                  pointsL0RadiusInf)
        
        self.LowestPosition = LowestPositionDegrees
        self.HighestPosition = HighestPositionDegrees
        return outline1, outline2, curvatureRadiusListUp, curvatureRadiusListLo

    def Conjugate_angular_circular(self, rb, rf, l1, l3sup, l3inf, beta):
        """It returns two list of points, the upper and the lower
        Keyword Arguments:
        rb    -- 
        l1    -- 
        l2sup -- 
        l2inf -- 
        l3sup -- 
        l3inf -- 
        beta  -- 
        """
        # curvature radius list
        curvatureRadiusPListUp = []
        curvatureRadiusPListLo = []
        curvatureRadiusListUp = []
        curvatureRadiusListLo = []
        # pressure angle
        pressureAngleListUp = []
        pressureAngleListLo = []

        offset1 = []
        offset2 = []

        # Making sure it works with radians
        LowestPositionDegrees = self.LowestPosition
        HighestPositionDegrees = self.HighestPosition
        self.LowestPosition = math.radians(self.LowestPosition)
        self.HighestPosition = math.radians(self.HighestPosition)
        
        # This list will contain all angles correspondeing to a peak
        peakRadiusOfCurvSup = []
        peakPointsSup = []
        # This list will contain all angles where there is no contact
        # with the follower
        noContactRadiusOfCSup = []
        noContactPointsSup = []
        # Angulos con rebaje
        negRadiusOfCurvSup = []
        negPointsSup = []
        # Angles with no good pressure angle
        badPressureAnglesSup = []
        badPressurePointsSup = []
        # All pressure angles that are problematic

        # Lower Follower
        # This list will contain all angles correspondeing to a peak
        peakRadiusOfCurvInf = []
        peakPointsInf = []
        # This list will contain all angles where there is no contact
        # with the follower
        noContactRadiusOfCInf = []
        noContactPointsInf = []
        # Angulos con rebaje
        negRadiusOfCurvInf = []
        negPointsInf = []
        # Angles with no good pressure angle
        badPressureAnglesInf = []
        badPressurePointsInf = []
        # All pressure angles that are problematic
        beta = math.radians(beta)

        # We suppose l2=l4=0
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        self.CreateYaxisAcceleration()
        outline1 = []
        outline2 = []
        # Fi1 is the angle that symbolize the rotation of the cam
        # Fi2 is the angle the follower does with the horizontal, Fi2
        # depends on Fi1
        numerator = math.pow(l1, 2) + math.pow(l3sup, 2) - math.pow(rb + rf, 2)
        denominator = 2 * l1 * l3sup
        # Fi2o is the initial angle the follower has because of the base
        # radius of the cam
        Fi2o = math.acos(numerator / denominator)
        for Fi1, ordinateAngle, angleDerivative, angleDerivative2 in zip(self.abcissaeValues,
                                                                         self.ordinatePosition,
                                                                         self.ordinateVelocity,
                                                                         self.ordinateAcceleration):
            angleDerivative = angleDerivative / self.AngularVelocity
            angleDerivative2 = angleDerivative2 / math.pow(self.AngularVelocity, 2)
            Fi2_1 = Fi2o + ordinateAngle
            Fi2_2 = (Fi2o - beta) + ordinateAngle
            cosFi2_1 = math.cos(Fi2_1)
            sinFi2_1 = math.sin(Fi2_1)
            cosFi2_2 = math.cos(Fi2_2)
            sinFi2_2 = math.sin(Fi2_2)
            # Vector OC calculation (the one that starts on the center of
            # the cam and ends on the center of the circular follower).
            # The initial basis is moving with the rotation
            # OC1_1 = (l1 + ((l2sup - l4sup) * sinFi2_1) - (l3sup * cosFi2_1))
            # OC2_1 = ((l2sup - l4sup) * cosFi2_1) + (l3sup * sinFi2_1)
            # OC1_2 = (l1 - ((l2inf + l4inf) * sinFi2_2) - (l3inf * cosFi2_2))
            # OC2_2 = ((-(l2inf + l4inf) * cosFi2_2) + (l3inf * sinFi2_2))
            OC1_1 = (l1 - (l3sup * cosFi2_1))
            OC2_1 = (l3sup * sinFi2_1)
            OC1_2 = (l1 - (l3inf * cosFi2_2))
            OC2_2 = (l3inf * sinFi2_2)

            Fi1 = math.radians(Fi1)
            # Now we change the basis of the vector
            cosFi1 = math.cos(Fi1)
            sinFi1 = math.sin(Fi1)
            OCx_1 = (cosFi1 * OC1_1) + (sinFi1 * OC2_1)
            OCy_1 = ((-sinFi1) * OC1_1) + (cosFi1 * OC2_1)
            OCx_2 = (cosFi1 * OC1_2) + (sinFi1 * OC2_2)
            OCy_2 = ((-sinFi1) * OC1_2) + (cosFi1 * OC2_2)

            offset1.append((OCx_1, OCy_1))
            offset2.append((OCx_2, OCy_2))
            
            # Tangentk vector calculation
            sinFi1PlusFi2_1 = math.sin(Fi1 + Fi2_1)
            cosFi1PlusFi2_1 = math.cos(Fi1 + Fi2_1)
            sinFi1PlusFi2_2 = math.sin(Fi1 + Fi2_2)
            cosFi1PlusFi2_2 = math.cos(Fi1 + Fi2_2)
            
            # vectorTanx = (((-l1) * sinFi1) +
            #               (((l2sup - l4sup) * cosFi1PlusFi2_1) *
            #               (1 + angleDerivative)) +
            #               (l3sup * sinFi1PlusFi2_1 * (1 + angleDerivative)))
            # vectorTany = (((-l1) * cosFi1) -
            #               ((l2sup - l4sup) * sinFi1PlusFi2_1 *
            #                (1 + angleDerivative)) +
            #               (l3sup * cosFi1PlusFi2_1 * (1 + angleDerivative)))
            vectorTanx = (((-l1) * sinFi1) +
                          (l3sup * sinFi1PlusFi2_1 * (1 + angleDerivative)))
            vectorTany = (((-l1) * cosFi1) +
                          (l3sup * cosFi1PlusFi2_1 * (1 + angleDerivative)))
            vectorTanx_2 = (((-l1) * sinFi1) +
                            (l3sup * sinFi1PlusFi2_2 * (1 + angleDerivative)))
            vectorTany_2 = (((-l1) * cosFi1) +
                            (l3sup * cosFi1PlusFi2_2 * (1 + angleDerivative)))
            
            # We now have to rotate the tangent vector -90 degrees to obtain
            # the normal vector. The rotation matrix for a counterclockwise
            # spin is: [S] = [[0, 1];[-1, 0]]
            vectorNormalx = vectorTany
            vectorNormaly = -vectorTanx
            vectorNormalx_2 = vectorTany_2
            vectorNormaly_2 = -vectorTanx_2
            vectorNormalUnitx = (vectorNormalx /
                                 math.hypot(vectorNormalx, vectorNormaly))
            vectorNormalUnity = (vectorNormaly /
                                 math.hypot(vectorNormalx, vectorNormaly))
            vectorNormalUnitx_2 = (vectorNormalx_2 /
                                   math.hypot(vectorNormalx, vectorNormaly))
            vectorNormalUnity_2 = (vectorNormaly_2 /
                                   math.hypot(vectorNormalx, vectorNormaly))
            # we now calculate OJ, the vector that points
            # to every point that generates the outline
            OJx_1 = (OCx_1 + (rf * vectorNormalUnitx))
            OJy_1 = (OCy_1 + (rf * vectorNormalUnity))
            OJx_2 = (OCx_2 + (rf * vectorNormalUnitx_2))
            OJy_2 = (OCy_2 + (rf * vectorNormalUnity_2))
            OJ1 = (OJx_1, OJy_1)
            OJ2 = (OJx_2, OJy_2)
            outline1.append(OJ1)
            outline2.append(OJ2)

            # Cheking the radius of curvature of the upper follower
            C = (-l1 + (l3sup * math.cos(Fi2_1) * (1 + angleDerivative)))
            D = l3sup * math.sin(Fi2_1) * (1 + angleDerivative)
            Cp2 = math.pow(C, 2)
            Dp2 = math.pow(D, 2)
            curvatureRadiusP = (math.pow((Cp2 + Dp2), (3 / 2)) /  # numerator
                                ((((C * l3inf *
                                    angleDerivative * cosFi2_2) + Dp2) *
                                  (1 + angleDerivative)) +
                                 Cp2 -
                                 (l1 * sinFi2_2 * l3sup * angleDerivative2)))
            curvatureRadius = curvatureRadiusP - rf
            if (curvatureRadius < 0):
                if abs(curvatureRadius) <= rf:
                    noContactRadiusOfCSup.append((curvatureRadius,
                                                  curvatureRadiusP))
                    noContactPointsSup.append(OJ1)
                else:
                    pass
            elif (curvatureRadiusP > 0):
                if (curvatureRadiusP == rf):  # There is a peak
                    peakRadiusOfCurvSup.append((curvatureRadius,
                                                curvatureRadiusP))
                    peakPointsSup.append(OJ1)
                elif (curvatureRadiusP < rf):  # Rebaje
                    negRadiusOfCurvSup.append((curvatureRadius,
                                               curvatureRadiusP))
                    negPointsSup.append(OJ1)
                else:
                    pass
            curvatureRadiusPListUp.append(curvatureRadiusP)
            curvatureRadiusListUp.append(curvatureRadius)
            
            # Lower follower
            CInf = (-l1 + (l3sup * math.cos(Fi2_2) * (1 + angleDerivative)))
            DInf = l3sup * math.sin(Fi2_2) * (1 + angleDerivative)
            Cp2Inf = math.pow(CInf, 2)
            Dp2Inf = math.pow(DInf, 2)
            curvatureRadiusPInf = (math.pow((Cp2Inf + Dp2Inf), (3 / 2)) /  # numerator
                                   ((((CInf * l3inf *
                                       angleDerivative * cosFi2_2) + Dp2Inf) *
                                     (1 + angleDerivative)) +
                                    Cp2Inf -
                                    (l1 * sinFi2_2 * l3inf * angleDerivative2)))
            curvatureRadiusInf = curvatureRadiusPInf - rf
            # print(curvatureRadiusInf)
            if (curvatureRadiusInf <= 0):
                if abs(curvatureRadiusInf) <= rf:
                    noContactRadiusOfCInf.append((curvatureRadiusInf,
                                                  curvatureRadiusPInf))
                    noContactPointsInf.append(OJ2)
                else:
                    pass
            elif (curvatureRadiusPInf > 0):
                if (curvatureRadiusPInf == rf):  # There is a peak
                    peakRadiusOfCurvInf.append((curvatureRadiusInf,
                                                curvatureRadiusPInf))
                    peakPointsInf.append(OJ1)
                elif (curvatureRadiusPInf < rf):  # Rebaje
                    negRadiusOfCurvInf.append((curvatureRadiusInf,
                                               curvatureRadiusPInf))
                    negPointsInf.append(OJ1)
                else:
                    pass

            curvatureRadiusPListLo.append(curvatureRadiusPInf)
            curvatureRadiusListLo.append(curvatureRadiusInf)

            # Time to calculate the pressure angle
            OCmodulusSup = math.hypot(OC1_1, OC2_1)
            numerator = (math.pow(OCmodulusSup, 2) +
                         math.pow(l3sup, 2) -
                         math.pow(l1, 2))
            denominator = 2 * OCmodulusSup * l3sup
            gamma = math.acos(numerator / denominator)
            numerator = (l3sup * math.sin(gamma) * angleDerivative)
            denominator = (OCmodulusSup -
                           (l3sup * math.cos(gamma) * angleDerivative))
            alpha = math.atan(numerator / denominator)
            pressureAngleSup = gamma + alpha - (math.pi / 2)
            if (math.degrees(pressureAngleSup) > 30):
                badPressureAnglesSup.append(math.degrees(pressureAngleSup))
                badPressurePointsSup.append(OJ1)
            else:
                pass

            pressureAngleListUp.append(math.degrees(pressureAngleSup))
            # Lower follower
            # Time to calculate the pressure angle
            OCmodulusInf = math.hypot(OC1_2, OC2_2)
            numerator = (math.pow(OCmodulusInf, 2) +
                         math.pow(l3inf, 2) -
                         math.pow(l1, 2))
            denominator = 2 * OCmodulusInf * l3inf
            gamma = math.acos(numerator / denominator)
            numerator = (l3inf * math.sin(gamma) * angleDerivative)
            denominator = (OCmodulusInf -
                           (l3inf * math.cos(gamma) * angleDerivative))
            alpha = math.atan(numerator / denominator)
            pressureAngleInf = gamma + alpha - (math.pi / 2)
            if (math.degrees(pressureAngleInf) > 30):
                badPressureAnglesInf.append(math.degrees(pressureAngleInf))
                badPressurePointsInf.append(OJ2)
            else:
                pass

            pressureAngleListLo.append(math.degrees(pressureAngleInf))
        # Falta el tratamiento final
        ShowCurvatureRadiusFailuresConjugateCircular(peakRadiusOfCurvSup,
                                                     negRadiusOfCurvSup,
                                                     noContactRadiusOfCSup,
                                                     peakPointsSup,
                                                     negPointsSup,
                                                     noContactPointsSup,
                                                     peakRadiusOfCurvInf,
                                                     negRadiusOfCurvInf,
                                                     noContactRadiusOfCInf,
                                                     peakPointsInf,
                                                     negPointsInf,
                                                     noContactPointsInf)

        ShowPressureAngleFailuresConjugateCircular(badPressureAnglesSup,
                                                   badPressurePointsSup,
                                                   badPressureAnglesInf,
                                                   badPressurePointsInf)
        self.LowestPosition = LowestPositionDegrees
        self.HighestPosition = HighestPositionDegrees

        return (outline1, outline2, curvatureRadiusPListUp,
                curvatureRadiusPListLo, curvatureRadiusListUp,
                curvatureRadiusListLo, pressureAngleListUp,
                pressureAngleListLo, offset1, offset2)

    def Desmo_plain(self, rb, betha, excentricity, separation):
        """Returns the outline as a list of points of a demodromic cam
        
        arguments:
        separation: A float that represents the distance between followers

        excentricity: A float that represents the excentricity
        """
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        sinBeta = math.sin(betha)
        cosBeta = math.cos(betha)
        initDistance = (rb + (excentricity * sinBeta)) / cosBeta
        outlineCurrent = []
        for fi1, position, deriv in zip(self.abcissaeValues,
                                        self.ordinatePosition,
                                        self.ordinateVelocity):
            cosFI = math.cos(math.radians(fi1))
            sinFI = math.sin(math.radians(fi1))
            # cosFIplusPI = math.cos(math.radians(fi1) + math.pi)
            # sinFIplusPI = math.sin(math.radians(fi1) + math.pi)
            # aplying change of basis
            x = (((deriv * cosBeta) * cosFI) +
                 ((((position + initDistance) * cosBeta) -
                   (excentricity * sinBeta)) * sinFI))
            y = ((((deriv * cosBeta) * (-sinFI)) +
                 ((((position + initDistance) * cosBeta) -
                   (excentricity * sinBeta)) * cosFI)))
            cosFIpi = math.cos(math.radians(fi1) + math.pi)
            sinFIpi = math.sin(math.radians(fi1) + math.pi)
            point1 = (x, y)  # dudas acerca del signo de la derivada
            x2 = (((-deriv * cosBeta) * cosFIpi) +
                  ((((-position - initDistance + separation) * cosBeta) -
                    (excentricity * sinBeta)) * sinFIpi))
            y2 = ((((-deriv * cosBeta) * (-sinFIpi)) +
                   ((((-position - initDistance + separation) * cosBeta) -
                     (excentricity * sinBeta)) * cosFIpi)))
            point2 = (x2, y2)
            outlineCurrent.append(point1)
            outlineAnti.append(point2)

        return outlineCurrent, outlineAnti, separation

    def Desmo_circular(self, rb, excentricity, separation, followerRadius):
        """Returns a the curve and the anti
        Keyword Arguments:
        rb             -- 
        excentricity   -- 
        separation     -- 
        followerRadius -- 
        """
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        initDistance = math.sqrt(math.pow((rb + followerRadius), 2) -
                                 math.pow(excentricity, 2))
        pointList = []
        pointListA = []
        for angle, position, velocity in zip(self.abcissaeValues,
                                             self.ordinatePosition,
                                             self.ordinateVelocity):
            # Whe have a rotating basis so we multiply the two components
            # of the vector with the change of basis matrix
            # [S] = [[cosFI, sinFI];[-sinFI, cosFI]]
            # the rotation is clockwise
            cosFI = math.cos(math.radians(angle))
            sinFI = math.sin(math.radians(angle))
            OCx = (cosFI * excentricity) + (sinFI * (position + initDistance))
            OCy = ((-sinFI) * excentricity) + (cosFI *
                                               (position + initDistance))
            # OC is the point that goes from the center of the cam to the
            # center of the follower's wheel
            # We calculate tangent vector
            vectorTanx = (((velocity - excentricity) * sinFI) +
                          ((position + initDistance) * cosFI))
            vectorTany = (((velocity - excentricity) * cosFI) -
                          ((position + initDistance) * sinFI))
            # We now have to rotate the tangent vector -90 degrees to obtain
            # the normal vector. The rotation matrix for a counterclockwise
            # spin is: [S] = [[0, 1];[-1, 0]]
            vectorNormalx = vectorTany
            vectorNormaly = -vectorTanx
            # We have interest in the unit vector
            vectorNormalUnitx = (vectorNormalx /
                                 math.hypot(vectorNormalx, vectorNormaly))
            vectorNormalUnity = (vectorNormaly /
                                 math.hypot(vectorNormalx, vectorNormaly))
            # we now calculate OJ, the vector that points
            # to every point that generates the outline
            OJx = (OCx + followerRadius * vectorNormalUnitx)
            OJy = (OCy + followerRadius * vectorNormalUnity)
            OJ = (OJx, OJy)
            pointList.append(OJ)
            
            # The other half of the problem
            cosFIpi = math.cos(math.radians(angle) + math.pi)
            sinFIpi = math.sin(math.radians(angle) + math.pi)
            OCxA = ((cosFIpi * excentricity) +
                    (sinFIpi * (-position - initDistance + separation)))
            OCyA = (((-sinFIpi) * excentricity) +
                    (cosFIpi *
                     (-position - initDistance + separation)))
            
            vectorTanxA = (((-velocity - excentricity) * sinFIpi) +
                           ((-position - initDistance + separation) * cosFIpi))
            vectorTanyA = (((-velocity - excentricity) * cosFIpi) -
                           ((-position - initDistance + separation) * sinFIpi))
            # Time to rotate the vector
            vectorNormalxA = vectorTanyA
            vectorNormalyA = -vectorTanxA
            
            # We have interest in the unit vector
            vectorNormalUnitxA = (vectorNormalxA /
                                  math.hypot(vectorNormalxA, vectorNormalyA))
            vectorNormalUnityA = (vectorNormalyA /
                                  math.hypot(vectorNormalxA, vectorNormalyA))
            # we now calculate OJ, the vector that points
            # to every point that generates the outline
            OJxA = (OCxA + followerRadius * vectorNormalUnitxA)
            OJyA = (OCyA + followerRadius * vectorNormalUnityA)
            OJA = (OJxA, OJyA)
            pointListA.append(OJA)

        return pointList, pointListA

    def Desmo_angular_plain(self, rb, l1, l2sup, l2inf):
        """
        Keyword Arguments:
        rb -- 
        """
        outline = []
        dsideOne = math.radians(self.ordinatePosition[0])
        dsideTwo = math.asin(-math.sin(dsideOne) + ((l2inf - l2sup) / l1))
        self.delta = math.degrees(dsideOne - dsideTwo + math.pi)
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        Fi2o = math.degrees(math.asin((rb - l2) / l1))
        for Fi1, ordinateAngle, rawDerivative in zip(self.abcissaeValues,
                                                     self.ordinatePosition,
                                                     self.ordinateVelocity):
            ordinateDerivative = ((rawDerivative / self.AngularVelocity) *
                                  math.radians(self.beta) / self.beta)
            Fi2 = Fi2o + ordinateAngle
            sinFi2 = math.sin(math.radians(Fi2))
            cosFi2 = math.cos(math.radians(Fi2))
            l3 = (l1 * cosFi2) / (1 + ordinateDerivative)
            OJ1 = (l1 + (l2 * sinFi2) - (l3 * cosFi2))
            OJ2 = ((l2 * cosFi2) + (l3 * sinFi2))
            # Time to change the basis to the fix one
            sinFi1 = math.sin(math.radians(Fi1))
            cosFi1 = math.cos(math.radians(Fi1))
            OJx = ((OJ1 * cosFi1) + (OJ2 * sinFi1))
            OJy = ((OJ1 * (-sinFi1)) + (OJ2 * cosFi1))
            OJ = (OJx, OJy)
            outline.append(OJ)

        return outline

    def GenerateCircularAnti(self, rb, excentricity, followerRadius,
                             separation):
        """Returns a the curve and the anti
        Keyword Arguments:
        rb             -- 
        excentricity   -- 
        separation     -- 
        followerRadius -- 
        """
        # Quitar el radio del seguidor
        # followerRadius = 0
        # excentricity = 0
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        initDistance = math.sqrt(math.pow((rb + followerRadius), 2) -
                                 math.pow(excentricity, 2))
        pointList = []
        pointListA = []
        for angle, position, velocity in zip(self.abcissaeValues,
                                             self.ordinatePosition,
                                             self.ordinateVelocity):
            # The other half of the problem
            cosFIpi = math.cos(math.radians(angle) + math.pi)
            sinFIpi = math.sin(math.radians(angle) + math.pi)
            OCxA = ((cosFIpi * excentricity) +
                    (sinFIpi * (-position - initDistance + separation)))
            OCyA = (((-sinFIpi) * excentricity) +
                    (cosFIpi *
                     (-position - initDistance + separation)))
            
            vectorTanxA = (((-velocity - excentricity) * sinFIpi) +
                           ((-position - initDistance + separation) * cosFIpi))
            vectorTanyA = (((-velocity - excentricity) * cosFIpi) -
                           ((-position - initDistance + separation) * sinFIpi))
            # Time to rotate the vector
            vectorNormalxA = vectorTanyA
            vectorNormalyA = -vectorTanxA
            
            # We have interest in the unit vector
            vectorNormalUnitxA = (vectorNormalxA /
                                  math.hypot(vectorNormalxA, vectorNormalyA))
            vectorNormalUnityA = (vectorNormalyA /
                                  math.hypot(vectorNormalxA, vectorNormalyA))
            # we now calculate OJ, the vector that points
            # to every point that generates the outline
            OJxA = (OCxA + followerRadius * vectorNormalUnitxA)
            OJyA = (OCyA + followerRadius * vectorNormalUnityA)
            OJA = (OJxA, OJyA)
            pointListA.append(OJA)

        return pointListA

    def GenerateOscilatingPlaneAnti(self, rb, l1, l2sup, l2inf, c):
        """Generates the outline for the cam with an angular follower
        with a plane instead of a point as the contact part
        Keyword Arguments:
        rb -- 
        l1 -- 
        l2 --
        """
        # This list will contain all angles whose radius value is zero
        anglesE0Radius = []
        pointsE0Radius = []
        # This list will contain all angles whose radius value
        # is less than zero
        anglesL0Radius = []
        pointsL0Radius = []

        LowestPositionDegrees = self.LowestPosition
        HighestPositionDegrees = self.HighestPosition
        self.LowestPosition = math.radians(self.LowestPosition)
        self.HighestPosition = math.radians(self.HighestPosition)
        # Passing values to radiants
        outline = []
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        self.CreateYaxisAcceleration()
        Fi2o = math.asin((rb - l2sup) / l1)
        Fi2o = 0
        for Fi1, ordinateAngle, ordinateDeriv, ordinateDeriv2 in zip(self.abcissaeValues,
                                                                     self.ordinatePosition,
                                                                     self.ordinateVelocity,
                                                                     self.ordinateAcceleration):
            Fi2 = Fi2o #+ ordinateAngle
            Fi2inf = math.asin(c - math.sin(Fi2))
            #print(math.degrees(Fi2inf))
            delta = Fi2 - Fi2inf + math.pi
            print(math.degrees(delta) + Fi1)
            firstNumerator = (math.pow(math.cos(Fi2), 2) *
                              (c - math.sin(Fi2)) *
                              math.pow(ordinateDeriv, 2))
            firstDenominator = math.pow((1 -
                                         math.pow(
                                             c - math.sin(Fi2), 2)), 3 / 2)
            firstAddend = firstNumerator / firstDenominator
            secondNumerator = ((math.sin(Fi2) * math.pow(ordinateDeriv, 2)) -
                               (math.cos(Fi2) * ordinateDeriv2))
            secondDenominator = math.sqrt((1 -
                                           math.pow(
                                               c - math.sin(Fi2), 2)))
            secondAddend = secondNumerator / secondDenominator
            # ordinateDeriv2 = firstAddend + secondAddend
            # ordinateDeriv = -((math.cos(Fi2) /
            #                   math.sqrt(1 - math.pow((c -
            #                                           math.sin(Fi2)), 2))) *
            #                   ordinateDeriv)
            sinFi2inf = math.sin(Fi2inf)
            cosFi2inf = math.cos(Fi2inf)
            # from here
            l3 = (l1 * cosFi2inf) / (1 + ordinateDeriv)
            OJ1 = (l1 + (l2sup * sinFi2inf) - (l3 * cosFi2inf))
            OJ2 = ((l2sup * cosFi2inf) + (l3 * sinFi2inf))
            # Time to change the basis to the fix one
            sinFi1 = math.sin(math.radians(Fi1) + delta)
            cosFi1 = math.cos(math.radians(Fi1) + delta)
            OJx = ((OJ1 * cosFi1) + (OJ2 * sinFi1))
            OJy = ((OJ1 * (-sinFi1)) + (OJ2 * cosFi1))
            OJ = (OJx, OJy)
            outline.append(OJ)
            # Careful with the derivatives ................
            # For this cam, only is necessary to confirm
            # if the Radius of curvature is correct
            firstAddend = ((1 + (2 * ordinateDeriv)) /
                           math.pow((1 + ordinateDeriv), 2)) * sinFi2inf
            secondAddend = (ordinateDeriv2 /
                            math.pow((1 + ordinateDeriv), 3)) * cosFi2inf
            giantExpression = firstAddend + secondAddend
            radiusOfCurvature = l2inf + (l1 * giantExpression)
            if (radiusOfCurvature == 0):  # the outline has a peak
                anglesE0Radius.append(radiusOfCurvature)
                pointsE0Radius.append(OJ)
            elif (radiusOfCurvature < 0):
                anglesL0Radius.append(radiusOfCurvature)
                pointsL0Radius.append(OJ)
            else:
                pass
        # Return the possible faliures
        ShowCurvatureRadiusFailures(anglesE0Radius,
                                    anglesL0Radius,
                                    pointsE0Radius,
                                    pointsL0Radius)
            
        self.LowestPosition = LowestPositionDegrees
        self.HighestPosition = HighestPositionDegrees
        return outline

    def GeneratePlain(self, rb, beta, excentricity, separation):
        """Returns the anti part of a constant width desmodromic cam
        
        arguments:
        separation: A float that represents the distance between followers

        excentricity: A float that represents the excentricity
        """
        # This list will contain all angles whose radius value is zero
        anglesE0Radius = []
        pointsE0Radius = []
        # This list will contain all angles whose radius value
        # is less than zero
        anglesL0Radius = []
        pointsL0Radius = []
        # Not finished yet more testing remains
        self.CreateYaxisPosition()
        self.CreateYaxisVelocity()
        sinBeta = math.sin(math.radians(beta))
        cosBeta = math.cos(math.radians(beta))
        #sinBeta = 0
        #cosBeta = 1
        initDistance = (rb + (excentricity * sinBeta)) / cosBeta
        outlineAnti = []
        for fi1, position, deriv, deriv2 in zip(self.abcissaeValues,
                                                self.ordinatePosition,
                                                self.ordinateVelocity,
                                                self.ordinateAcceleration):
            s = position + initDistance
            seb = (s - (excentricity * (sinBeta / cosBeta)))
            sebAnti = separation - seb
            sebDeriv = -deriv * cosBeta
            cosFIpi = math.cos(math.radians(fi1) + math.pi)
            sinFIpi = math.sin(math.radians(fi1) + math.pi)
            x2 = ((sebDeriv * cosFIpi) +
                  (sebAnti * sinFIpi))
            y2 = ((sebDeriv * (-sinFIpi)) +
                  (sebAnti * cosFIpi))
            point = (x2, y2)
            outlineAnti.append(point)

            # checking if the outline is OK by checking the curvature radius
            curvatureRadius = ((position + initDistance + deriv2) *
                               cosBeta) - (excentricity * sinBeta)
            if (curvatureRadius == 0):  # the outline has a peak
                anglesE0Radius.append(curvatureRadius)
                pointsE0Radius.append(point)
            elif (curvatureRadius < 0):
                anglesL0Radius.append(curvatureRadius)
                pointsL0Radius.append(point)
            else:
                pass
        # Return the possible faliures
        ShowCurvatureRadiusFailures(anglesE0Radius,
                                    anglesL0Radius,
                                    pointsE0Radius,
                                    pointsL0Radius)


        return outlineAnti    

    def GeneratePlainAnti(self, rb, beta, excentricity, separation):
        """Returns the anti part of a constant width desmodromic cam
        
        arguments:
        separation: A float that represents the distance between followers

        excentricity: A float that represents the excentricity
        """
        # Curvature radius list
        curvatureRadiusList = []
        # This list will contain all angles whose radius value is zero
        anglesE0Radius = []
        pointsE0Radius = []
        # This list will contain all angles whose radius value
        # is less than zero
        anglesL0Radius = []
        pointsL0Radius = []
        # Not finished yet more testing remains
        # self.CreateYaxisPosition()
        # self.CreateYaxisVelocity()
        # self.CreateYaxisAcceleration()
        sinBeta = math.sin(math.radians(beta))
        cosBeta = math.cos(math.radians(beta))
        initDistance = (rb + (excentricity * sinBeta)) / cosBeta
        outlineAnti = []
        for fi1, position, deriv, deriv2 in zip(self.abcissaeValues,
                                                self.ordinatePosition,
                                                self.ordinateVelocity,
                                                self.ordinateAcceleration):
            deriv = deriv / self.AngularVelocity
            deriv2 = deriv2 / math.pow(self.AngularVelocity, 2)
            cosFIpi = math.cos(math.radians(fi1) + math.pi)
            sinFIpi = math.sin(math.radians(fi1) + math.pi)
            x2 = (((-deriv * cosBeta) * cosFIpi) +
                  ((((-position - initDistance + separation) * cosBeta) -
                    (excentricity * sinBeta)) * sinFIpi))
            y2 = ((((-deriv * cosBeta) * (-sinFIpi)) +
                   ((((-position - initDistance + separation) * cosBeta) -
                     (excentricity * sinBeta)) * cosFIpi)))
            point = (x2, y2)
            outlineAnti.append(point)

            # checking if the outline is OK by checking the curvature radius
            newS = separation - position - initDistance
            curvatureRadius = ((newS - deriv2) *
                               cosBeta) - (excentricity * sinBeta)
            if (curvatureRadius == 0):  # the outline has a peak
                anglesE0Radius.append(curvatureRadius)
                pointsE0Radius.append(point)
            elif (curvatureRadius < 0):
                anglesL0Radius.append(curvatureRadius)
                pointsL0Radius.append(point)
            else:
                pass
            curvatureRadiusList.append(curvatureRadius)
        # Return the possible faliures
        ShowCurvatureRadiusFailures(anglesE0Radius,
                                    anglesL0Radius,
                                    pointsE0Radius,
                                    pointsL0Radius)

        return outlineAnti, curvatureRadiusList
